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We develop a semi-classical approximation to electron spin resonance in quantum spin systems, 
based on the rotor or non-linear sigma model. The classical time evolution is studied using molec- 
ular dynamics while random initial conditions are sampled using classical Monte Carlo methods. 
Although the approximation may be especially powerful in two dimensions, we apply it here to one- 
dimensional systems of large spin at intermediate temperatures, in the presence of staggered and 
uniform magnetic fields. We first test the validity of the semi-classical approximation by comparing 
the magnetization to quantum Monte Carlo results on S — 2 chains. Then we calculate the ESR 
spectrum, finding broad coexisting paramagnetic and spin wave resonances. 

PACS numbers: 76.20. +q, 75.10.Pq, 03.65.Sq 
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I. INTRODUCTION 

Electron spin resonance probes the dynamics of inter- 
acting spin systems at zero wave-vector in a magnetic 
field. The intensity of adsorption of microwave radia- 
tion of frequency w is a LuG^: a (uj) where G„ a (w) is the 
retarded Green's function of = J2 r ^r> the a compo- 
nent of the total spin operator. Since [Sj,, Hq] — for an 
SU(2) invariant Hamiltonian, Ho, such as the Heisenberg 
model, the ESR spectrum remains exactly the same as 
that for a single spin. Thus ESR is a highly sensitive 
probe to small anisotropies in the Hamiltonian. In order 
to discuss the effect of small anisotropies, the calculation 
of ESR spectra has to be very precise. An approxima- 
tion, that may be reasonable to study other properties 
such as inelastic neutron scattering spectrum, might lead 
to an errorneous change of ESR, even in the absense of 
anisotropy. Then it could not be applied to ESR. This is 
a challenging aspect of ESR theory, especially in strongly 
correlated systems. 

Older theories of ESR 1 -^ are based on high tempera- 
ture expansions at T > J, the exchange coupling, where 
each spin behaves approximately independently or else 
on spin- wave theory at T < Tjy, the Neel temperature, 
typically of order J. For the two dimensional Heisen- 
berg model on the square lattice, Tjv = 0, leaving a large 
temperature range where neither of these approaches to 
ESR applies. In this "renormalized classical region" an 
approach to dynamics was developed by Chakravarty, 
Halperin and Nelson^ (CHN) based on the two dimen- 
sional classical rotor model (CRM). The CRM is defined 
with two vector variables rii and Li at each site, and the 
Hamiltonian 



constraints 



H 



lD-2 



2x«_ 



where D is the dimension of the lattice (D = 2 in Ref. |3J) 
and b is the lattice constant, rii and Lj are subject to 



m = 1, 
rii ■ Li = 0. 



(2) 
(3) 



Physically, Li represents the angular momentum of the 
rotor rii, and the two parameters p s and Xu± represent, 
respectively, the spin stiffness and the transverse compo- 
nent of uniform susceptibility. 

The dynamics of the system is defined by the Hamil- 
tonian ([T]) together with the Poisson brackets 



7 



(4) 



K>f} = o, 



where a, (3, 7 = 1, 2, 3 and e Q/37 is the Levi-Civita symbol. 

The classical time evolution is calculated using molec- 
ular dynamics simulations starting from random ini- 
tial conditions which are generated with classical Monte 
Carlo methods. CHN and Tyc, Chakravarty and 
Halperin 4 applied this method to calculate the neutron 
scattering cross-section at wave- vector near the antifcrro- 
magnetic point [tt/ol, n/a), where a is the lattice constant 
of the spin system. The purpose of this paper is to ex- 
plore the applicability of this method to dynamics at zero 
wave-vector, which is relevant to ESR. 

In fact, as shown in Refs.d and@, this classical method 
can also be applied in dimension D — 1, at least for suffi- 
ciently large spin magnitude, S. While Neel order doesn't 
occur even at T — for the one dimensional Heisenberg 
model there is another important characteristic low en- 
ergy scale, A of order 



Je 



-ITS 



(5) 



For integer S this is the energy of the gapped triplet 
magnons, predicted by Haldane and following from the 
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behavior of the quantum 0(3) non- linear sigma model 
(NLSM), basically the continuum limit of Eq. (1). At 
low temperatures, T«A ESR can be calculated^ from 
transitions between single magnon states (or, if appropri- 
ate symmetry breaking is included in H, from magnon 
production processes). Magnon-magnon interactions be- 
come unimportant in this temperature range because the 
magnons are dilute, with density oc e~ A ^ T . Instead, 
at higher temperature, the magnons are dense, and the 
NLSM behaves rather classically. The rotor model is 
nothing but a lattice version of the NLSM; thus CRM 
describes the classical dynamics of the NLSM. 

The present approach to ESR based on the CRM can 
be used in the intermediate temperature range 

A<T<J5 2 , (6) 

where the magnons are dense and the CRM is still valid 
as an effective theory. Unfortunately, this temperature 
range may not exist for the spin S = 1 chain, for which 
the Haldane gap is known to be A ~ 0.41J. However, the 
gap becomes smaller for higher 5 as in eq. ([5]) and the 
temperature range ([6]) becomes well-defined for higher 
spins, perhaps starting at 5 = 2, for which the Haldane 
gap is already as small 9 as A~ 0.089 J. 

While the classical method has a wider range of T and 
5 over which it is applicable in D = 2, we focus on the 
D = 1 case in this paper, due to the computational cost of 
the classical technique, which appears to be quite severe 
for ESR studies. We hope to return to the D = 2 case 
in the future. Therefore, in this paper, we discuss ESR 
in ID systems in the previously unexplored temperature 
regime ©. We will discuss conditions for the present 
approach to be applicable, in more detail, in Sec. MI A[ 

Following Refs. |j and la, initial states are generated 
by classical Monte Carlo simulation of the 0(3) CRM at 
given temperature. The real-time correlation function is 
then obtained by solving the classical equation of mo- 
tion for the 0(3) CRM with the initial condition. ESR 
spectrum can be obtained from Fourier transform of the 
real-time correlation function. 

Among various possible anisotropics, in this paper we 
discuss the staggered field, which is known to have most 
interesting effects on ESR spectra in 5 = 1/2 and 5 = 1 
cases. The effective Hamiltonian including the staggered 
field is given by 

H = J2 Js i ■ s i+i - HS j - (-lYhSf, (7) 

3 

where H and h is the uniform field and the transverse 
staggered field respectively. Throughout this paper, we 
set gfis = ft — ks = 1 for simplicity. 

While the staggered field might seem unphysical, it 
often does appear effectively in actual quantum antifcr- 
romagnetic chains. When an external magnetic field is 
applied to a material with a staggered crystal structure 
along the chain, the staggered field is effectively gener- 
ated through a staggered g tenso r 10 ' 11 and also through 



a staggered Dzyaloshinskii-Moriya interaction. 12 Exam- 
ples of such materials include the typical 5=1 Hal- 
dane chain material Ni(C 2 H 8 N2)2N0 2 C104 (NENP)±2, 
and the 5 = 1/2 Heisenberg antiferromagnetic chain Cu 
benzoateii. 

For ESR in the 5 = 1/2 antiferromagnetic chain at 
low temperatures, field theory approach leads to diverg- 
ing linewidth of the paramagnetic peak at lower tem- 
peratures, and an appearance of a new peak when the 
temperature is lowered further down to zero 15 ' 16 . For 
ESR in the 5 = 1 Haldane chain at low temperatures, 
the violation of selection rule by the staggered field leads 
to the appearance of a new peak at the frequency equal 
to the Haldane gap, corresponding to creation of a sin- 
gle magno n 10 ' 11 . These approaches are only justified in 
a "quantum" regime and no longer valid for a Haldane 
chain at the intermediate temperatures which is the 
focus of the present paper. 

The 5 = 2 Heisenberg antiferromagnetic chain, for 
which the present approach would be relevant, is not just 
a theoretical toy model. Granroth et at 17 reported an ex- 
perimental evidence of Haldane gap in an 5 = 2 antiferro- 
magnetic Heisenberg chain compound MnCls(bpy) (bpy 
= bipyridine). MnCL^bpy) is similar to NENP. It has a 
quite small single ion anisotropy D/J < 0.04 ± 0.02 and 
has a staggered crystal structure. The staggered crys- 
tal structure would cause the staggered g tensor, which 
produces an effective staggered field when an external 
magnetic field is applied. Thus it would be interesting to 
measure ESR in MnCl3 (bpy) and compare to the present 
theory. 

Our results may be applied to a wider range of one- 
dimensional systems such as quantum spin ladders as 
well, since they can also be described by the same CRM 6 . 
For spin ladders, there are various possible generaliza- 
tions of the staggered field. The results of the present 
paper can be directly applied when the staggered field is 
unfrustrate d 18 ' 19 , for example as in 

n 

« = E E [ JS ^ • s^+i - ss*j (-ir +j hsi P ] 

/*=1 3 
n—1 

+ E E J±s nj ■ s h+ij> ( 8 ) 

M=l 3 

where n is the number of legs and /i is the leg index. In 
this case, the staggered field can be handled in the same 
way as in the chain (J7J. 

In fact, even for gapless one-dimensional systems such 
as half-integer spin Heisenberg antiferromagnetic chains, 
A defined in Eq. ([5]) is a characteristic energy scale. At 
energy scales of order A, the system renormalizes from 
the weak coupling regime, where classical methods can 
be used, to the non-trivial critical point induced by the 
topological term in the effective Lagrangian. At temper- 
atures small compared to A the low temperature theory 
for the 5=1/2 Heisenberg mode l 15 ' 16 can be used. (For 
the 5=1/2 case, this theory is valid at any tempera- 
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ture T <C J.) On the other hand, in the temperature 
range ([5]), the present approach based on the CRM is 
valid even in gapless systems. 

Physically, ESR in the presence of the staggered field 
provides an interesting case of crossover of dynamics be- 
tween two different regimes. When the effect of the stag- 
gered field is weak, the ESR spectrum is dominated by 
the paramagnetic resonance at uj ~ H, and the staggered 
field causes its broadening and shift. However, when the 
effect of the staggered field is strong, the system is or- 
dered along the staggered field. The ESR spectrum is 
dominated by the spin-wave type fluctuation around the 
ordered state. This is similar to the case in which the 
system has a Neel order spontaneously, but is different 
in that the "order" is imposed externally by the stag- 
gered field. Nevertheless, theory of antiferromagnetic 
resonanc o 20 i 21 developed for the spontaneously ordered 
state can be modified and applied to the present case, as 
discussed for example in Ref. [l0|. It describes the limit of 
the strong staggered field, where the imposed Neel order 
is perfect. 

In the case of a spontaneous ordering, a phase tran- 
sition separates the ordered and disordered phases. In 
the present case, there is no phase transition but only a 
smooth crossover between the two regimes. The descrip- 
tion of the crossover is generally much more difficult than 
that of the limiting cases. 

For the 5 = 1/2 chain, the field theory approach suc- 
cessfully gives the broadening and shift of the paramag- 
netic resonance, when the staggered field is small (or the 
temperature is sufficiently high) . It also explains the ESR 
spectrum in the low temperature limit, when the system 
is ordered along the staggered field. The antiferromag- 
netic resonance in this case is renormalized due to strong 
quantum fluctuations in the 5 = 1/2 chain. This can 
be well described in terms of elementary excitations of 
quantum sine-Gordon field theory. However, despite the 
integrability of the quantum sine-Gordon field theory, dy- 
namical susceptibility at finite temperature has not been 
obtained exactly. Thus the theoretical description of the 
crossover between two regimes still remains unsolved, al- 
though a numerical result based on exact diagonalization 
was reported 2 ^. 

In other systems, the crossover is even less understood. 
The present approach can, where it is valid, numerically 
describe the nontrivial crossover of the ESR spectrum 
between two regimes as we will demonstrate. 

This paper is organized as follows. In Sec. HH we intro- 
duce the quantum 0(3) NLSM as an effective low-energy 
theory for the Hamiltonian ([7]), and review its renormal- 
ization. ESR spectrum is also formulated in terms of the 
0(3) NLSM. We will then discuss the classical approxi- 
mation of the 0(3) NLSM, including the range of validity 
of the approximation, in Sec. lIIII The CRM is introduced 
as a lattice version, suitable for numerical calculation, of 
the classical limit of the 0(3) NLSM. In Sec. |TV]we will 
obtain ESR spectra by numerically solving the equation 
of motion for the CRM. We will develop a theory of anti- 



ferromagnetic spin- wave resonance in Sec. |V1 which can 
be identified with the new resonance at higher frequency 
found in the numerical results. Finally in Sec. I VII we 
summarize the paper. 



II. 0(3) NONLINEAR SIGMA MODEL 
A. Definition of the model 

Let us introduce the 0(3) NLSM as an effective field 
theory of the antiferromagnetic Heisenberg chain, and 
summarize its properties. Although we will eventually 
treat the system classically, first we need to clarify the 
effects of quantum fluctuation in order to determine the 
appropriate parameters in the effective classical model. 

The 0(3) NLSM is defined in terms of two fields 
n(x) and L(x) which are related to the original spin Sj 
through 



Sj = (-iySn(x)^l-(^P^ + aL{x). (9) 

Here x — ja, and a is a lattice spacing. Hereafter we 
set the lattice spacing to a = 1. n(x) and L(x) satisfy 
constraints 

n 2 {x) = 1, (10) 
n(x) ■ L(x) = 0, (11) 

which are necessary to keep the constraint S 2 = 5(5 + 1). 
Then n(x) and L(x) satisfy the following commutation 
relations 

[L a (x),L^(y)} =ie a ^L^(x)S(x-y), (12) 
[L a (x), n?(y)} = ie a ^(x)S(x - y), (13) 
[n a (x),nP(y)}=0. (14) 

0(3) NLSM is derived after substituting j| to and 
taking the continuum limit. The Hamiltonian density of 
0(3) NLSM is given by 

^=^~ (5,n) 2 + f L 2 -HL z -^n^ (15) 
2g 2 gc 

where c = 2J5 is the spin- wave velocity, g = 2/5 is the 
coupling constant, and A^ is defined as 

A h = VUSh, (16) 

Later in this Section, we will show that A/,, is the 
staggered-field-induced gap. The coupling constant g is 
actually subject to renormalization owing to quantum 
fluctuations, which we will discuss in the next subsec- 
tion, g = 2/5 should be understood as the bare coupling 
constant. 
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Here we introduce useful normalization which measure 
the energy scale in units of c = 2 JS. 



H' = 

r = 

H' = 

h' = 

K = 
t' = 



u 

2JS 
T 

2JS 
H 

2JS 
h 

2JS 
A h 

2JS 
2JSt 

U) 

2JS 



= V2h> 



(17a) 

(17b) 

(17c) 

(17d) 

(17e) 
(17f) 
(17g) 



We specify the normalized, dimensionless parameters by 
adding the prime symbol as in the above. In this nota- 
tion, the Hamiltonian (|15[) is rewritten as 



1 



2ff 



9 -L 2 - H'L Z - ^ 
2 9 



-^-n x . (18) 



Because Eq. (IT51) (or Eq. (|18p) is harmonic in L(x), 
L(x) — H/gc = n x (d t n + H x n)/gc = n x (d t >n + 
H' x ri)/g can be integrated out. Then the Lagrangian 
is given by 

jSf' = d t ,n ■ (d t >n + H' x n) — Jt' 



H' 



^(d.n) 2 + ^n x , (19) 
25 2# 



with the constraint (j 10|) . In this expression, it is clear 
that the coupling constant g represents the degree of 
quantum fluctuation, since it plays the role of Planck's 
constant. 

Although there is no energy scale in the Lagrangian 
density, 0(3) NLSM at H = h = has an excitation 
gap, which corresponds to the Haldane gap A. In field 
theory language, the mass gap is dynamically generated. 
The low-energy phenomena of the 0(3) NLSM can be 
effectively described by the field theory of triplet bosons 
with the mass A: 



Jzf 



2.9 



1 
2g~c 



A 2 n 2 (x), 



(20) 



where A^ = J H' 2 6 a , z + A' h . The field m = (0,n«,n z ) 
contains the components of n transverse to the staggered 
field direction. We assumed n x (t, x) w 1 because of large 
h. m represents the elementary excitation in the large h 

regime and has a dispersion E'j} a ^ = \jk 2 + A^ 2 where 
k is the wave number and a — y, z. This dispersion 
relation indicates that the gap for the z component is 
different from the gap for the y component: 



K = 



H' 2 



(22) 
(23) 



where A^ is the staggered-field- induced gap (117e[) . The 
gap (|22p is derived within the classical approximation ig- 
noring quantum fluctuation. On the other hand, in the 
absence of the staggered field, the quantum fluctuation 
leads to the Haldane gap A. Roughly speaking, the ex- 
citation gap of the system is given as 



A' ~ max(A', A' h ) 



(24) 



In ESR measurements, we apply an oscillating field 
with frequency uj and polarization a. In this paper we 
consider Faraday configuration where the polarization is 
perpendicular to the direction of the uniform field. ESR 
absorption intensity 



(25) 



is written in terms of the dynamical susceptibility 
Xaa(w). According to the linear response theory, 
the imaginary part of Xaa(w) is related to Fourier 
component of the retarded Green function G^ a (t) = 
-W(i)([i¥(t),i?(0)]>, namely X 'M - -ImG£». 
The retarded Green function G^ a (uj) is defined as 



1 1 dte tut {[L%(t),L%(0)}}. 
la 



(26) 



= J dxL a (x) is the total magnetization. In general, 
the ESR spectrum thus depends on the polarization a. 

In Ref. [la, it was shown that for the system ([7]) with 
the staggered field in x-direction, the polarization de- 
pendence can be exactly determined by the equation of 
motion: 



H 2 cos 2 $ + w 2 sin 2 $ „ 



X;'», (27) 



where $ is the angle between the polarization a and x- 



without the constraint (|10|). Hereafter (d^n) means 
(d^n) 2 = (l/c 2 )(d t n) 2 - (d x n) 2 = {d v nf - (d x n) 2 . 

When the applied field H becomes larger than A while 
the staggered field h = 0, the gap is closed. On the other 
hand, when staggered field h is large, the excitation gap 
is still open even if H > A because the Lagrangian (|T9]) 
is approximated as follows. 



Jzf' 



1 

2^ 



A' 2 

y K) 2 



2.9 



z\2 



2.9 



(m z ) 



(21) 



Furthermore, ESR spectrum for circular polarization 
can be related to 

/>oo 

x';_H=Rc/ dte^([L+(t),L T (0)}) 7 (28) 
Jo 

where L ± = L x ± iL y . The discussion in Ref. |l6| can be 
extended to circular polarization to obtain 

x'|-M = + ( 29 ) 
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That is, ESR spectra for different polarization in Fara- 
day configuration are related to each other, and thus it 
is sufficient to compute the spectrum for one particular 
case. 

In this paper, we consider the circular polarization 
case, eq. ([28jl . A short calculation yields 



(w) = (f-e-"/ T )Re / dte iut {L+{t)L^(0)), (30) 



which is convenient for our purpose since eq. ((30 j) is well- 
defined even if L a (t) is a classical vector. In this paper, 
we numerically solve the classical equation of motion to 
obtain the dynamical correlation function 



(£+(t)2^(Q)). 



(31) 



This gives the desired ESR spectrum by eq. (f3"U)) . Our 
numerical approach to this problem will be discussed in 
Sec. Eland Appendix [5] 



Using Eq. (33) for g R (E), Eq. (06) becomes: 



d In h' D 



N-3 



1 



dlnE 2(N - 2) \n(E / A) 
This gives the renormalized staggered field 



h' R (E) = h' 



ln(£/A) 
In (J/ A) 



(JV-3)/2(W-2) 



(38) 



(39) 



Namely, for the case of our interest TV = 3, the renormal- 
ization of the staggered field h! is absent. The renormal- 
ized Lagrangian is given by just replacing g by g R (E) in 
eq. (TT9)) . Because h' R (E) = h' is kept unchanged (within 
the leading order considered in this paper) for N = 3, 
there is no logarithmic correction in the staggered-ficld- 
induced gap (ITB1) 



III. 



CLASSICAL APPROXIMATION 



B. Renormalization group transformation 



A. Conditions for the classical approximation 



Let us consider the renormalization group transforma- 
tion of general O(N) NLSM whose Lagrangian is given 
by 

7r(dt,n + H' xn) 2 -^(d x n) 2 + — n x . (32) 

This bare theory is defined at the energy scale J. we 
renormalize this theory down to an energy scale E < J. 

When H = h = 0, the effective coupling gii(E) at 
energy scale E satisfies 2 ^ 



f 



9 9R{E) 



N-2, (E 

In — 

2tt \ J 



(33) 



As the bare coupling satisfies 1/g = (N — 2) ln(J/A)/27r, 
we can simplify the renormalized coupling constant 



1 



N-2 E 

— - In — . 

2tt A 



9r(E) 

For the case of our interest, N = 3, we find 

, , 2?r 
9r{E) 



ln(£/A)' 



(34) 



(35) 



Next we consider the renormalization of the staggered 
field. We ignore the influence of h! on the renormalization 
of <?, but consider the renormalization of h' produced by 
<7, at lowest order. This is determined by the anomalous 
dimension of the field r» 23 i 24 . 



dh\ 



R 



-l{g R {E))h' R 



d\nE 

To one loop order, the anomalous dimension 7(g) is 

N-3 



7(5) 



-in 



-9 



(36) 



(37) 



In the Introduction, we have briefly discussed the con- 
dition (J6j) for the classical approximation being justified. 
Here we discuss the condition in some more detail, with 
several additional conditions. 

The latter condition in eq. ©, 



T< JS 2 



(40) 



is required for the validity of the effective theory in the 
continuum. It is equivalent to the requirement that the 
correlation length of the antiferromagnetic order param- 
eter is much longer than the lattice spacing. In fact, as 
pointed out in Ref. [|| two different temperature regimes 
can be distinguished for large S. 



and 



T < ^max ~ 2 JS, 



(41) 



(42) 



In the regime the 0(3) NLSM description is valid. 
On the other hand, when eq. (|4"2"j) holds, the quantum 
spin chain may be directly approximated as a classical 
spin chain and then the continuum description is ap- 
plied. The resulting dynamics is equivalent to that of 
the classical 0(3) NLSM, although the effective param- 
eters are estimated in a different manner. For simplic- 
ity, throughout the rest of this paper, we assume the 
regime (|4II) . Although ESR in the higher temperature 



range T, 



(2) 



JS 2 < T is also an interesting nontrivial 



problem, it is out of scope of the present paper. 

ESR is studied under an applied magnetic field H . If 
the ground state is fully polarized along H, 0(3) NLSM 
is not appropriate as an effective theory. Thus, our ap- 
proach also requires that the uniform field should be 
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much weaker than its saturation field, which is order of 
JS: 



H < JS. 



(43) 



The first inequality in eq. © is necessary for the classi- 
cal approximation to hold. Namely, that the temperature 
higher than the gap, implies high density of thermally ex- 
cited magnons. This leads to a breakdown of the quan- 
tum mechanical picture^ of ESR based on transitions of 
independent magnons. In fact, as we have discussed in 
Sec. Ill Al the staggered field induces a gap as in eq. ([T^|. 
Since actual gap of the system is given as eq. (IM1) . the 
first inequality in eq. (|6|) should be replaced by 



A = max(A, A h ) < T. 



(44) 



Namely, the temperature must be higher than not only 
the Haldane gap A, but also the staggered-field induced 
gap Ah. Together with eqs. (|4"0")) and (fToT) . this requires 



T 2 JS 3 

h<< TTs << —- 



(45) 



Later, we will discuss an antifcrromagnetic spin- wave the- 
ory of ESR. It is justified when 



TT 2 T 2 



h > 



US(ln(T/A)Y 



(46) 



where we have assumed eq. (|4"5j) . For a larger spin S, 
the Haldane gap A is expected to become exponentially 
small as in eq. ([5]). Thus the range of h satisfying both 
eqs. (|4"5"|) and (|46l) becomes wide for large S. In fact, we 
will demonstrate that the spin-wave theory prediction 
agrees very well for a S = 10 chain with staggered field, 
in a range of parameters. On the other hand, for 5 = 2, 
when the classical approximation is valid, the spin-wave 
theory is not quite justified. Correspondingly, a broad 
peak is observed instead of a sharp resonance. Neverthe- 
less, the broad peak may be understood as a remnant of 
the spin- wave resonance. 

For discussion of ESR based on the classical dynamics 
of the NLSM, the frequency uj of the applied oscillating 
field also should be much lower than the temperature. 

u < T (47) 

This also implies 

H < T, (48) 
because ESR absorption usually occurs for oj > H . 

B. Effective Hamiltonian 

In the Introduction, we argued that the system should 
be described by the classical 0(3) NLSM in the interme- 
diate temperature range ([5]). 



Even in the temperature range (|6]), we cannot simply 
ignore quantum fluctuations at energy scale above the 
temperature T. Their effects are taken into account by 
the rcnormalization group. The effective classical Hamil- 
tonian for the 0(3) NLSM may be obtained by using the 
renormalized parameters, eqs. (|35|) and (f39j) . and setting 
the energy scale to the temperature: E = T. 

On the other hand, we do not consider the renormaliza- 
tion of H ■ L term since J dx L is a conserved quantity 
when h — 0. Therefore any renormalization that does 
occur should vanish at h = and would be negligible at 
small h. The renormalized Hamiltonian density is thus 



given as 



jfa = —{d x nf + ^L 2 - HL Z 



c 
29~R 



2cg R 



^-n x , (49) 



where g R — g R (T). In terms of dimensionless parame- 
ters, it reads 

= ^-(d x nf + ^-L 2 H'L* —n\ (50) 
25fl 2 g R 

As we will discuss below, the coefficient of the (d x n) 2 
term and the L 2 term in the classical Hamiltonian (|4"9")l 
can be identified respectively with £T/2 and 1/2x„_l, 
where £ is the classical correlation length of n and 
Xu± = (3/2)xu is proportional to the zero-field uniform 
susceptibility \u, both at zero fields H = h = 0. Namely, 



e= 2^T ln A 
1 T 

X "^2^ ln A 



(51) 
(52) 



Thus, the effective Hamiltonian density (|49|) can be also 
written as 

■aSi = ^(d x n) 2 + -^—L 2 (x) - HL z (x) - X u±A 2 h n x , 
2 2x„_l 

(53) 

which is the form used in Ref. [6J. 

In numerical calculations we discretize the effective 
field theory (IBTfl) and consider the Hamiltonian, 



1 



b 2 g R 



■j+i - 



— Al 



9R 

(54) 

on a lattice with the lattice spacing b. This is nothing 
but the CRM (TTJ in dimension D = 1. We note that it is 
not necessary to take b equal to the lattice spacing a = 1 
of the original spin chain; the CRM may be regarded as 
a lattice regularization of the classical, continuum 0(3) 
NLSM. Usually b > a is taken, because the eq. (|54]l is in- 
troduced to describe long-distance assymptotic behavior 
of the spin system. The system size L = j\fb is propor- 
tional to the number of rotors Af. We use Af = 16 for 
our numerical calculations. As we will see in FIG. [TJ 
TV = 16 is large enough to reproduce consistent values 



7 



of the magnetization density MJL with quantum Monte 
Carlo calculations and low-field expansion. This is be- 
cause the correlation length of rij is much shorter than 
L due to the relatively high temperature T ~ J. 

Let us now demonstrate that the above identifications 
are valid within the classical theory. To do so, we assume 
the Hamiltonian density in the form of eq. (|53l) , and then 
show that £ and Xu are indeed the correlation length and 
the uniform susceptibility. To calculate the correlation 
function of n, it is convenient to integrate out L and 
obtain the Hamiltonian in terms of n. For H — h = 
Ah = 0, it reads 



(55) 



The last expression is the NLSM Hamiltonian on a dis- 
cretized one-dimensional lattice with the lattice constant 
a. In fact, it is equivalent to the classical Heisenberg 
chain. The equilibrium statistical properties of the clas- 
sical Heisenberg chain are studied by Fisher—. The cor- 
relation function (n a (x)n b (0)) of the classical Hamilto- 
nian ([55)l is obtained exactly by using transfer matrix 



method, (n a (x)n b (0)} = (l/3)6 ab exp(-\x\/£), whi ch 



dicates that Eq. (|51[) is indeed the correlation length of 
the order parameter n. 

Next we introduce an infinitesimal uniform field H. 
The zero-field uniform susceptibility Xu is defined as 
Xu = hrtiH-j-o dM/dH where M = J dx(L z ) is the total 
magnetization. We calculate the uniform susceptibility 
Xu of the classical NLSM ([55]). with h = A h = 0. After 
integrating out L, the Hamiltonian becomes 



H= dx 



T^fdn 
~2~ \dx 



(n J 



(56) 



Here we defined n 1 - = (n x ,n v ,0). The uniform suscepti- 
bility is thus obtained as 

Xu = M / dxdn^) 2 ) = r,Xu±, (57) 

where () is the thermodynamic expectation value at H = 
0. (n x ) 2 = 2/3 follows from the isotropy of the Hamilto- 
nian and the constraint (| 10|) . Thus, (|52|) is proved to be 
a transverse component of the uniform susceptibility Xu- 



C. Magnetization 

Our primary interest in this paper is in the dynamics 
of the system. However, before going into the dynamics, 
it would be important to establish the validity of the 
present approach by considering static properties. This 
was done earlier in Ref. [2(| in the zero field limit. 

Here we demonstrate its validity for nonzero magnetic 
field H, which is relevant for ESR, by discussing the mag- 
netization and its dependence on H and T. In a classical 



system with U(l) symmetry around H (thus with h = 0), 
we find the following interesting identity 

M=^ J dxdy{L*(x)L x (y)) ch (58) 

where M = (i^)ci is the classical uniform magnetiza- 
tion. This identity is valid for any H and T. The proof 
of ([58|) is given in the Appendix [A] We can easily con- 
firm that this classical magnetization is approximately 
equal to the uniform magnetization of the corresponding 
quantum chain. 

Since we approximate the renormalization group equa- 
tion at the lowest order, the uniform field H and the 
staggered field h do not affect the renormalization of the 
coupling constant g, except that H gives the energy scale 
E for the renormalization if H > T . At small but fi- 
nite field, the uniform susceptibility has a small deviation 
from Q57j) ■ which is proportional to H 2 . We expand the 
uniform magnetization M with respect to small H : 



AI 

T 



Xu±H 



Xu±H 2 



dx((n^(x)) 2 (n^(0)) 2 ) Q 



(-)o denotes the thermal average by the classical zero-field 
Hamiltonian (151)1) . The second term can be calculated by 
using the following formula for four-point function, 



(n a {x)n b (x)n c {Q)n d {Q)), 



1 P -3M/£ 

_ S ab S cd + (3(6 ac S bd + S ad S bc ) - 26 ab 5 cd ), 

9 45 v ' 



and the result is, 



M 

T 



= Xu±H 



4 tXu±H 2 
135 T 



0{H 5 



(59) 



The uniform (differential) susceptibility is given by 



:XuA 



2 Xu^H 2 



3-7TC 



In 



15 T 

max (T, H) 
A 



cH 2 
15ttT 2 



In 



(60) 

2 / max (T, H) 



A 



(61) 



We employ the energy scale E = max (T, H) in the log 
correction in (|59l) and (|6"Tj) instead of E = T as in the 
previous Section, because the magnetic field can exceed 
the temperature when we do not consider ESR. If H > T, 
we should take a cut-off scale as E — H. 

We have calculated the magnetization density in the 
quantum S — 2 antiferromagnetic chain by quantum 
Monte Carlo simulation, and also in the effective CRM by 
classical Monte Carlo (CMC) simulation. Our method of 
CMC simulations of the 0(3) CRM is explained in Ap- 
pendix [B] The quantum Monte Carlo simulation was 
done using the codes provided by ALPS projec t 27 ' 28 . In 
the 5 = 2 chain, the magnetization density is defined 
simply by M/L = (S^)/L. The numerical results are 
compared with the low-field expansion ([5T)f in Fig. [TJ 
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The magnetization density in the quantum S = 2 chain 
agrees well with that in the effective classical 0(3) NLSM, 
and with eq. (f59|) . We note that there is no adjustable 
parameter in this comparison, owing to the fact that the 
magnetization is a conserved quantity. For higher field 
H, there is a visible discrepancy between the numerical 
results and the analytical prediction (1591) . This is pre- 
sumably due to the higher order terms in H which are 
ignored in eq. ([59)) , and not because of breakdown of the 
classical description. In fact, the CMC result for the ef- 
fective classical 0(3) NLSM agrees quite well with the 
quantum 5 = 2 chain, even at higher field H . 

In Fig. [5] we show the difference between the magne- 
tization density calculated by the CMC simulations and 
the low field expansion (|50[) . The difference between the 
two results is indeed proportional to H 5 , which is the 
next order in the expansion. These results on the magne- 
tization supports the validity of the classical description 
for S > 2 chain, also in a finite magnetic field H . 
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FIG. 2. The open circles are the classical Monte Carlo data 
(T/J — 1) which are shown in Fig. [T] The filled circles are 
the difference of CMC data and (|59[1 . This subtracted data 
agree well with the solid line 0.0022 x H whose coefficient is 
determined by fitting. 



D. Equations of motion 
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1.5 



0.5 



CMC (T/J =1, N=16) 
CMC (T/J =0.5, N=16) 
QMC (T/J =1, N=20) 
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1st order 
3rd order (T/J =1) 
3rd order (T/J =0.5) 




1 1.5 2 

Magnetic field H/J 

FIG. 1. The rescaled magnetization density M/L = M/L x 
l/(JX«i.) f° r 0(3) NLSM is plotted. The open circles and tri- 
angles are numerical data by classical Monte Carlo simulation. 
The filled circles and triangles with error bars are obtained by 
quantum Monte Carlo simulation. We found finite-site effects 
in QMC data very small, by comparing simulations on 20 sites 
and 40 sites. We set the unit of the vertical axis so that the 
simple relation M'/L w (2/3) H' holds at low H' = H/J < 1. 
The solid line represents the low field expansion of M/L up 
to 1st order of H, which is M'/L = M/(LJ X u±) = \(H / J) 
by (|59p . The dashed and dotted lines are 3rd order approx- 
imations of (|59|) at temperature T/J — 1 and T/J = 0.5 
respectively. The low H expansion (|59jl is consistent with 
the numerical data in the regime H/T < 1. We emphasize 
that CMC and QMC data are consistent in any value of the 
magnetic field. 



The classical dynamics of the CRM (|54| is studied by 
solving the classical equation of motion. In the classical 
theory, the commutation relations (|T2l , (fl3l) and (| L4[) are 
replaced by the Poisson brackets ((4|). These brackets lead 
to the equations of motion, 



drij 

dL 3 

dt' 



(g R (T)L,-H')x njl 
1 



rij x (n J+ i + n,-_i) 



H' X La 



a; 2 

9r(T)' 



(62) 



(63) 



IV. ESR SPECTRUM 

In the present approach, the ESR spectrum is obtained 
from the classical dynamics of the effective 0(3) NLSM 
theory. The classical dynamical correlation function is 
calculated as follows. First we generate initial states us- 
ing the classical Monte Carlo Method so that the proba- 
bility distribution of initial states is identical to the Boltz- 
mann weight with the Hamiltonian (|53[) at given temper- 
ature. For each initial state, we solve the equations of 
motion (1621) and (|63[) numerically. We must pay careful 
attention to the total energy, which is a conserved quan- 
tity. Determination of the ESR lineshape requires the 
asymptotic, long-time behavior of the dynamical corre- 
lation function. The time evolution was obtained up to 
time t' max = 2000 with the time step 5t' = 0.001. 

Naive numerical integration of the equations of mo- 
tion results in violation of the energy conservation. This 
makes the scheme unsuitable for ESR calculation, for 
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which, as we have discussed in the Introduction, high ac- 
curacy is required. We applied symplectic method a 29 ' 30 
to assure the conservation law. We will give detailed 
explanations about these numerical methods in Ap- 
pendix [Bj The initial state is generated by classical 
Monte Carlo simulation of the 0(3) NLSM at the given 
temperature. 




2000 



FIG. 3. Time evolution of L+(t)L^(0) for S = 2, T' = 
0.15, and H' = 0.005, obtained by solving classical equation 
of motion, for one initial state generated by the classical 
Monte Carlo simulation. The time evolution appears chaotic 
for nonzero staggered field h', before taking ensemble average. 



An example of the time evolution of L^(t)L^(0) ob- 
tained by numerically solving the equation of motion, 
for one initial state is shown in Fig. [3] For zero stag- 
gered held, its exact solution is given by harmonic os- 
cillation. However, in the presence of a nonvanishing 
staggered field, the time evolution looks chaotic. After 
Fourier transform, the spectrum is also noisy reflecting 
the "noise" in the time evolution. 
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FIG. 4. The Fourier transform \(LfL T }\(u}') of the dynamical 
correlation function pT)l for S = 2, T' = 0.15, H' = 0.005, 
and h' — 0.008, obtained by averaging over various number 
of sample initial states. The approximately white noise is 
reduced by increasing the number of samples. 



The correlation function (|3Tj) is defined by thermal en- 
semble average over initial states. In our calculation, 
the above steps are repeated for 10 4 Monte Carlo sam- 
ples of initial states, and the average is taken. After 
taking the average, we obtain the dynamical correlation 
(L^(i)L^ (0)). The averaging over many samples of the 
initial states reduces the noise in the spectrum. In Fig. 0] 
we show the absolute value of the dynamical correlation 
function (|3T1) at h' = 0.008, obtained by average over 
1, 100, 1000, and 10000 samples of initial states. We note 
that the data contains both signal and noise. For 10000 
samples, the result shows small fluctuation due to noise, 
around the signal. For 1 sample, on the other hand, the 
observed data is dominated by the noise, which randomly 
takes positive or negative values. (In Fig. 01 the absolute 
value is shown.) 

In the range of frequency < u>' < 0.15 we are inter- 
ested, the power spectrum of the noise is approximately 
independent of the frequency u>'; it may be regarded as 
a white noise. The reduction of the noise by averag- 
ing is clear in the figure, and the noise is proportional 
to 1/ y/Ng, where N s is the number of initial state sam- 
ples. However, for a realistic number of samples (we used 
N s ~ 10 4 ) the noise is still not completely negligible. 
This effectively constrains resolution in the frequency, as 
we discuss below. 

Because of the finite time interval t' max , resolution in 
the frequency space is at most 2ir/t' mlix ~ 0.003. It is 
well known that a simple Fourier transform of the time 
dependence up to the cutoff time t' max leads to an artifi- 
cial spreading of resonances known as spectral leakage. 31 
To suppress the spectral leakage, a window function is 
multiplied to the data before Fourier transform. Here we 
apply the Gauss window function! 32 ' 33 Namely, instead 
of (L+{t)L^(0)), we take the Fourier transformation of 



(64) 



(L+(t)L T (0))exp(- — 



Obviously, the width of the Gaussian window a must be 
smaller than t' m&K . Usually the width a is still taken as 
the same order as t' max , for example a ~ QAt' max . How- 
ever, we find that the obtained spectrum is affected by 
the noise for such a choice of a. In order to reduce the 
noise, we take a much smaller than t' max . The reduction 
of the noise by windowing is discussed in Appendix [C] 
There is a trade-off between the reduction of noise (better 
for larger rr) and resolution in the frequency space (bet- 
ter form smaller rj). In this paper, we choose a = 100, 
which corresponds to the resolution rx^ 1 w 0.01 in the 
frequency space. 

ESR spectra obtained numerically are shown in 
Fig. ©for S = 10, and in Fig. © for S = 2. They are 
related to the spectrum discussed in Fig. [4] via eqs. [251 
and [301 We note that at lu 1 ~ is suppressed by 

the factor a/ 2 compared to \(Lj,L^)\(cj'). 

The dependence on the spin quantum number S comes 
only through the effective coupling constant gn, as given 
by eqs. ([331) an d ©• From Fig. [5] and Fig.© we can see 
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that resonance splits to two peaks as h is increased. One 
peak at L) rs H results from the paramagnetic resonance. 
The intensity of this absorption peak becomes smaller as 
the staggered field becomes larger. For S — 10, the sec- 
ond peak at higher frequency is sharp, while the original 
paramagnetic peak at uj' s» H' is almost invisible except 
for h' = 0.5 x 10~ 3 . On the other hand, for S = 2, the 
original paramagnetic peak persists and the second peak 
is very broad for the studied parameter range. 

The second peak at higher frequency lo > H must be 
caused by the staggered magnetic field, because it be- 
comes dominant while the paramagnetic resonance peak 
is suppressed, as we increase the staggered field. As we 
increase S, this peak survives while the paramagnetic 
resonance vanishes as shown in Fig. [5j In the next sec- 
tion, we will discuss the physical origin of the new peak 
at u > H, as well as the difference between S = 10 
(Fig. [5]) and S = 2 (Fig. The shape of each peak is 
asymmetric, suggesting that the lineshape is not a sim- 
ple Lorcntzian, in contrast to many cases with strong 
exchange interactions! 1 ' 2 ' 15 

We note that the location of the paramagnetic peak in 
Fig.[B]is at u/ w 0.015 which is higher than expected H' = 
0.005. This is an artifact due to the windowing (|64l) . as 
discussed in Appendix [Q 
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FIG. 5. ESR spectra for the S = 10 case. We use dimen- 
sionless parameters T' = 0.3, H' = 0.15. While the paramag- 
netic resonance peak vanishes for h! > 1.5 x 10~ 2 , the peak 
in high frequency side survives. The high-frequency peak is 
much sharper compared to the corresponding resonance in the 
S = 2 case (Fig.©. 
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FIG. 6. ESR spectra for the dimensionless staggered field h' = 
h/2JS = 0, 5.0 x 10~ 3 , 8.0 x 10~ 3 and 1.1 x 10" 2 are shown. 
Here we assume S — 2. The dimensionless uniform field and 
temperature are respectively given as H' — H/2JS = 0.005 
and T' — T/2JS = 0.15. In addition to the paramagnetic 
resonance (h! = 0) at cj ~ H , we can observe the very broad 
peak in higher frequency side. 



served in Figs. [S]and[n] ■ Here we discuss linearized fluc- 
tuations around the antiferromagnetic order externally 
imposed by the staggered field. In this aspect, it is dis- 
tinguished from the standard theory of antiferromagnetic 
resonance, in which the antiferromagnetic order is caused 
by a spontaneous symmetry breaking. 

We initiate our discussion by taking the Lagrangian in 
the large h regime, which is (|2Tj) . because the additional 
peak originate in the large h limit h — > +oo. By replacing 
the bare parameters to the renormalized ones, we obtain 
the classical spin wave theory: 



2<m(T) 



{d t ,mf + (d x m) 2 



+ A; 2 (m y ) 2 + Af(m z ) i 



(65) 



This Hamiltonian indicates that the classical spin dynam- 
ics in the large h limit is governed by the two harmonic 
modes m v and m z . These oscillating modes have the 
eigenfrequencies (|2"2l and (1231) . 

Now let us discuss the condition for the spin- wave the- 
ory to be justified. The spin-wave theory is based on 
the assumption that the field n is polarized along the 
staggered field direction, and the fluctuation around the 
polarized groundstate is small. Therefore, the condition 
can be written as 



V. SPIN WAVE THEORY 

We propose an antiferromagnetic spin wave theory in 
order to explain the second peak at higher frequency ob- 



(m v ) 2 + {m z ) 2 < 1. 



(66) 



Since A' z > for any h', we can expect ((m z ) ) < 
((m y ) 2 ). Thus, it is sufficient to require ((m v ) 2 ) <C 1. 
The Gaussian Hamiltonian (ltI5|) leads 



({myf)»g R {T)r 



2a: 



< i 
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and this immediately results in 



9RV T ) 



(67) 



In terms of the physical staggered field h, this condition 
can be written as eq. (|46[) . 

We note that, the classical approach requires ((6|), 
namely 



r'>max(A , ,A' h ,i/')- 



(68) 



In particular, T" A^ is required. In the classical 
limit 5* — > oo, the renormalized coupling constant gFt{T) 
approaches zero. This indicates that the temperature 
range, where both (|57)l and (j6"5[) hold, becomes larger as 
we increase the spin quantum number S. 

In the classical picture, ESR corresponds to precession 
of magnetic moments. If there is no anisotropic interac- 
tion, the total magnetic moment £t precesses around the 
field H with the frequency lo = H without any dissipa- 
tion. The modes mentioned above seem to affect the dy- 
namics of Lt through the equation of motion (1521 . The 
additional peak reflects the cigcnfrcquency of m which is 
different from the paramagnetic resonance frequency H. 

In order to see this, the following identity on the dy- 
namical susceptibility \-\ ( w ) is useful: 



x+-M = 



2m 



([A,L?]) 



H (lo - H) 



+ XAU(w), (69) 



where A — A'^n^/gii(T). This formula is easily obtained 
by integrating the left hand side by parts. First and 
second terms in the right hand side have a singularity 
only at lu — H. The last term XAi A m (|69|) contributes 
to the additional singularity in X+-- 

The dynamical correlation function of n z can be eas- 
ily calculated within the effective spin-wave Hamilto- 
nian (I65[) . Thus the spin- wave approximation predicts 
an additional resonance at the frequency identical to 
eq. 



A' 



W 



A', 2 . 



(70) 



Let us compare this with the numerical results in Figs. [5] 
and[6j We note that, Figs. [5] and [6J are shown for to' < T', 
where the classical approximation would work. 

First we discuss the S — 10 case, which is highly classi- 
cal. The resonance peak observable in Fig.[5]for non-zero 
h! can be identified with the antiferromagnetic spin-wave 
resonance discussed above, since the original paramag- 
netic resonance peak disappears. The resonance frequen- 
cies for several values of h! are plotted in Fig. and 
compared with the theoretical prediction ([70]) . Both con- 
ditions ([6"7[) and are satisfied in the cases studied in 
Fig. [7] Thus the linear spin-wave approximation of the 
classical 0(3) NLSM should be valid. In fact, we find a 
very good agreement with the spin-wave theory predic- 
tion in this case. 



The disappearance of the paramagnetic peak and the 
sharp spin-wave resonance corresponding to eq. ([70)1 
may be understood as consequences of small fluctuation 
around the polarized state along the staggered field. For 
a large spin such as S = 10, the coupling constant g^. is 
small and thus eq. (|67p is easily satisfied, leading to small 
fluctuation (|66[) . Thus the spectrum is well described by 
the spin- wave theory which gives the sharp resonance at 
frequency (|70[) . The original paramagnetic resonance at 
lu ~ H corresponds to global precession of spins around 
the magnetic field H \\ z. When eq. (|57]) holds, the spins 
are polarized along the staggered field and thus the global 
precession cannot occur; the paramagnetic resonance is 
expected to vanish. This is consistent with the observed 
behavior in Fig. [3] 

Next we discuss the S = 2 case, in which the quan- 
tum fluctuations are stronger corresponding to the larger 
value of gn- In this case, as we have seen in Fig. [51 two 
peaks are observed. The lower frequency peak represents 
the paramagnetic resonance at lo ~ H, which is almost 
independent of the staggered field. As in the S = 10 
case, the higher frequency peak would be identified with 
the antiferromagnetic spin-wave. However, in this case, 
a quantitative analysis of the resonance frequency is not 
possible because the peak is very broad. The survival 
of the paramagnetic peak and the broadness of the anti- 
ferromagnetic spin-wave resonance for S = 2 are conse- 
quences of large fluctuation owing to large gu, in contrast 
to the S — 10 case discussed above. 

We note that, at a lower temperature, the fluctua- 
tion around the polarized state becomes smaller and the 
spin-wave theory holds better. However, for 5 = 2, the 
classical approximation is no longer valid in this regime. 
Nevertheless, the antiferromagnetic spin- wave resonance 
should exist also in the quantum regime; in fact, it exists 
even in the S = 1 chain with a staggered field at low 
temperatur e) 10 ' 1 1 . Our classical calculation describes the 
broadening of the antiferromagnetic resonance at higher 
temperatures. 

Finally, we note that, our theory predicts the spin- 
wave resonance frequency (|22p to be independent of tem- 
perature. This is also in agreement with our numeri- 
cal calculations (not shown). While the temperature- 
independence of the resonance frequency might seem ob- 
vious, it is owing to the lack of the renormalization of 
the staggered field for N = 3 as shown in eq. (|39|) . In a 
similar analysis based on the 0(N) NLSM with N > 4, 
the "spin-wave" resonance frequency should depend on 
the temperature through the logarithmic correction. 



VI. SUMMARY AND DISCUSSION 

We discussed ESR in the classical limit based on 
the 0(3) NLSM and the corresponding CRM. Our dis- 
cussion is valid in the classical temperature regime, 
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FIG. 7. Staggered field dependence of the resonance frequency 
in S = 10 system at T' = 0.3, H' = 0.15. The results agree 
quite well with the prediction of the spin wave theory (I70p 
(solid line). This implies that the resonance peak in the S = 
10 chain under a staggered field indeed corresponds to the 
antiferromagnetic spin wave. 



max(A, A h ,H) < T < JS 2 . Here the field theoretical 
picture by 0(3) NLSM holds, and the system behaves 
classically. However the microscopic parameters in bare 
0(3) NLSM are renormalized by quantum fluctuations. 
The dynamics of the 0(3) NLSM with renormalized pa- 
rameters was then analyzed in the classical limit. Actual 
numerical calculations were done for the lattice version 
of the NLSM, namely the CRM. 

We have demonstrated that numerically calculated 
ESR spectra show that the paramagnetic and antifer- 
romagnetic resonance peaks coexist in this intermediate 
temperature regime. The latter is characteristic of this 
regime since it disappears in both lower and higher tem- 
perature regions. We analytically showed that the an- 
tiferromagnetic resonance frequency is uj = \J H 2 
based on the linearized spin wave-theory (|65j). This 
agrees well with the numerical results for a large spin 
quantum number S = 10. For S = 2, owing to larger 
fluctuation, the spin-wave theory is not quite justified 
in classical regime. Nevertheless, a broad resonance is 
observed in the spectrum, which is identified with the 
(remnant of) antiferromagnetic spin-wave resonance. 

It is instructive to compare the present results with 
ESR in the 5=1/2 antiferromagnetic chain. Although 
the 5=1/2 chain at low temperature is described by a 
different theoretical approach (bosonization), as we have 
discussed in the Introduction, there is a qualitative sim- 
ilarity between the two cases. That is, the paramagnetic 
resonance at the frequency uj ~ H is broadened and even- 
tually disappears as the temperature is lowered. On the 
other hand, at lower temperatures, the new (antiferro- 
magnetic resonance) peak at higher frequency becomes 
dominant. 

There is an important difference in the frequency of the 



antiferromagnetic resonance. In both cases, it is given by 
eq. (|2"2"j) . However, the staggered- field- induced gap has 
different dependence on the staggered field h. In general, 
the staggered-field-induced gap is given by 



Ah oc h 2 - 



(71) 



where 7 is the anomalous dimension of the staggered 
magnetization. For the low-energy limit of the 5=1/2 
chain, the staggered magnetization has the anomalous 
dimension 1/2. Thus it follows 15 that 



Ah oc h 2 l\ 



(72) 



up to logarithmic corrections. In contrast, in the present 
case, the anomalous dimension is basically zero and thus 



A h oc h 1 / 2 , 



(73) 



which is consistent with cq. p6[) . This difference can be 
understood as an effect of strong quantum fluctuations 
in the 5=1/2 chain at low temperatures, on the anti- 
ferromagnetic resonanc o 20 i 21 . 

The lineshape of each peak in the present system ap- 
pears to be asymmetric and thus non-Lorentzian. This 
would also be a significant difference from the 5=1/2 
case, in which the broadening of the paramagnetic peak 
has Lorcntzian forniji 5 - 

In order to apply the present formulation to quantum 
spin systems, we need 5 > 2 in the case of a single chain. 
The temperature range max(A, A^, if) <C T <C JS 2 is 
not wide enough for the 5=1 chains. For the 5 = 2 
Heisenberg antiferromagnetic chain, we confirmed the 
validity of the classical 0(3) NLSM approach in a fi- 
nite magnetic field by calculating the magnetization den- 
sity with quantum Monte Carlo simulations. The results 
show a good agreement with the effective classical 0(3) 
NLSM. 

We thus expect that our approach is applicable to the 
5 = 2 Haldane chain compound MnCi3 (bpy) , which may 
bear an effective staggered field due to its staggered crys- 
tal structure. It would be interesting to measure ESR 
spectrum and compare to our prediction, especially the 
appearance of the antiferromagnetic resonance and its 
broadening. The present results can also be applied to 
spin ladders or tubes, if the system is described the the 
classical 0(3) NLSM, and if the staggered field is unfrus- 
trated as in eq. (J5J. 

Finally, let us comment on extension of the present ap- 
proach to ESR in 2 dimensional antiferromagnets. The 
present approach can be extended to 2 dimensions, as 
discussed in the Introduction. The temperature range 
for the classical approach is wider in 2 dimensions, as it 
is not limited from below by the Haldane gap. This is 
favorable for the classical approach in 2 dimensions. On 
the other hand, the increase of the computational cost 
in higher dimensions may be a problem for calculation 
of ESR spectra which requires high precision and res- 
olution. Since ESR in 2 dimensional systems is rather 
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little understood, it would be interesting to pursue this 
approach in 2 dimensions. We plan to attempt this as a 
next step. 
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Appendix A: Proof of (f58|) 



Let us define the partition function 



Z(H) = / @n@L*[[5{n{x) ■ L{x))S{n 2 (x) - l)e 



■n/T 



. ,:Ali 

The key observation is that, thanks to the rotation invari- 
ance of the Hamiltonian, the partition function depends 
just on the length of H, namely Z{H) = Z(\H\). Let us 
consider the two infinitesimal variations of H from Hq — 
(0, 0, H ). The first variation is H = (0, 0, H +SH). The 
partition function can be expanded in terms of SH as 



1 



1. Initial states at equilibrium 

The initial states are generated by the Monte Carlo 
method as follows. This method is also used to study 
static properties in Sec. IIII CI For the sake of simplic- 
ity, we will discuss in terms of continuum variables for 
simplicity. The actual calculation is done for the CRM 
on a lattice, for which (d x rij) 2 is to be replaced by 

First we eliminate L from the classical Hamiltonian 
(f54"|) . For simplicity, we denote n(t — 0,x) and L(t = 
0,x) as no(x) and Lq(x) respectively. The constraint 
no • Lq = decrease the degrees of freedom of L$. For 
instance, we eliminate Lg from the Hamiltonian, 



W — - 
cl " 2g R (T) 

1 



(X „ \2 ■ 9r{T) t2 , z 2b' 
{o x n ) H — L — H L n 7^n, 



-(d x n ) 2 + 



T , n o ■ L o 



H' 



2 

9R(T) 
2 

2b! 



no, 



° 9r{T) 

2 . ( n ' -^o) 2 



where Lq and are defined as Lq = (Lq,Lq,0) and 
n a = ( n 7 n o>0) respectively. Lq is distributed with 
Gaussian distribution with an average —Xu±HriQnQ, and 
the variance Txu± ■ We can integrate out Lq because the 



Z { H) = Z { Ho) [l + -^)o8H + 0((8Hf), (A2) ^^^^^1 



up to the first order in SH. Here (-)o means that the 
expectation value with H = Ho- The second variation 
we consider is H = (5H 1 ,0, Ho). Again the partition 
function can be expanded in terms of SH'. The first 
order term actually vanishes here because {L^)o = due 
to the symmetry. The expansion up to the second order 



I 



2s«(T) 



(d x no) 



H' 



{n z of 



2ti 



9R(T) 



n x (Bl) 



Z(H) = Z(H ) 1 



1 

2T 2 



mfMSH'f 



0((SH') 



We compare the variation of H by the two infinitesimal 
variation of if. \H\ w Hq + SH for the first variation and 
\H\ w Hq + {SH') 2 /2H Q for the second. The isotropy of 
the partition function ([AT) leads that SH = (SH') 2 /2H . 
Thus, we find the identity 



2 ^ 



(A4) 



which is valid for any temperature and for any magnetic 
field H . This result (fAH is nothing but ((55)) . 



Starting from an arbitrary configuration of {n(x)}, we 
thermalize it by classical Monte Carlo method based 
on the Wolff algorithm^. Lq is calculated by Lq = 
—H'nQno/gii{T) + R. R is a random vector whose dis- 
tribution is Gaussian with an average and a variance 
TxuJ_- The constraint n ■ L = determines (^q)- 



2. Calculation of dynamical correlation functions 

We need to solve the classical equations of motion (|52")) 
and (|63p numerically with the initial states produced in 
Appendix [Bl We should pay attention to the conserved 
quantities, the total energy and the uniform magnetiza- 
tion (if there is no staggered field) in our system. Or- 
dinary numerical methods for solving the equations mo- 
tion, for instance, Runge-Kutta method or the predictor- 
corrector method, result in the violation of the conser- 
vation law. This violation is caused by discretizations. 
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Even if the continuum equations of motion preserve the 
conservation laws, the discrctizcd version of the equations 
of motion do not necessarily preserve them. 

This problem is serious in our situation. ESR spec- 
trum requires asymptotic behavior of the dynamical cor- 
relation function, namely (L^(t)L^,(0)) at t — > +00. The 
configuration {n(t,x), L(t,x)} may be far from the real 
one which exactly conserves the energy once the conser- 
vation law is violated. 



a. Symplectic methods 



according to (|B4[) and BCH formula. (IB5|) is much more 
accurate than that of 4th order Runge-Kutta method. 
Let us consider the total energy at time t and denote it as 
E(t). First we prepare initial state and fix E(t = 0). By 
applying the symplectic method, E(t) — E(0) is order of 
(St) 5 , thus, E(t)-E(0) ~ (St) 5 . However, by the Runge- 
Kutta method, the error is increasing at least linearly 
with respect to t. Thus, E(t) — E(0) ~ t becomes quite 
large if we want to consider large t. 

Thus, if we apply the symplectic method, the total en- 
ergy conserves with high accuracy during the discretized 
time evolution. 



We overcame this difficulty by applying "symplectic" 
method a 29 i 30 . This algorithm is based on the Suzuki- 
Trotter decompositio n 35 ! 36 of integrators. Applications 
to classical spin systems are explained in detail in Refs.l29l 
andl3fj 

Let us consider the time evolution of a scalar f(t) for 
simplicity. We assume the equation of motion for f(t) as 



( 1 
dt 



F(f(t))=F 1 (f(t))+F 2 (f(t)). 



Further we assume that we can exactly solve the equa- 
tions 



4f 
dt 



(B2) 



for i = 1,2. Infinitesimal time evolution of / from t to 
t + St has an exponential form 

f(t + St) = e 5ti) f(t) = e st ^ 1+ ^f(t). (B3) 



Since we can solve (IB2|) exactly, we know the explicit 
form of f)i and f)2 and thus, e 5tt)1 and e 5tt)2 . The essence 
of the symplectic method is in approximating the opera- 
tor e St ^ 1+t)2 ^ by exactly known operators e Stt)1 and 1 



,5*1)2 



The fourth order approximation is known 



n 



e p i ,5tln/2 e p i <5i(, 2 



,Pi8tt> x /2 



(B4) 



where p\ = P2 = Pi = P5 = 1/(4 — 4 1 / 3 ) and p% — 
-4 1 / 3 /(4 — 4 1 / 3 ). We can extend the decomposition (|B4|) 
to 2nth order one with n > 3. By applying Baker- 
Campbell-Hausdorff (BCH) formula repeatedly, we can 
show that the right hand side in flH} is e Stt )+o(( St ) 6 ) . it 
is easily proved that 2n-th order decomposition of this 
type is (2n + l)-th order approximation. It follows im- 
mediately that the symplectic method is more accurate 
than other naive methods and that the conserved quanti- 
ties during the calculation do not deviate from the exact 
value. Let us discuss a specific example about the latter. 
In Hamilton systems, the operator f) is the Hamiltonian 
itself. Let us write the approximated f) in (|B4p as f)4. 
The difference f)4 — f) is bounded as 



f)4 — rj = const, x (5i) E 



(B5) 



b. Application to classical rotor model 

The application of this symplectic method is as follows. 
We decompose the equations of motion (|62j) and (f63|) : 



dfUj = gn(T)Lj x rij 
dt.Lj = 



H' X 71, 



dfLj = H' x rij 



(B6a) 



(B6b) 



dftij = 

(B6c) 



Each of these equations are exactly solvable. For in- 
stance, (|B6al) is solved as follows: 



n (t' + St') = ra^ (*0 cos (sii(T)<ft' \L 3 (t)\ 
L j (t' + St') = L l (t') 



We consider this step from t' to t' + St' as e St A . In the 
same way we define the time step e St B and e St c as 

H' x n (t') 

nj(t' + St') = nj (t') cos(HSt') H / y ' sm(HSt') 

L 3 (t' + St') = Lj(t') cos(HSt') - — ^J^Ll sm(HSt') 



and 

nj(t' + 5t') = nj(t') 



L J (t' + 5t') = L J (t') + St'(^—n J (t') 

\9r\J- ) 



2ti N 

x (nj+itf) +nj-i(t')) — e x x n^t') 

9RV 1 ) / 
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respectively. Actual time evolution is considered as 
e st (A+B+c) _ Second order decomposition of e 8t'(A+B+c) 

ig e 8t'{A+B+C) _ e St'A/2 e 8t'B/2 e 6t'C e St'B/2 e 8t'C/2_ J n tn j s 

manner, we can apply these to the symplectic method. 
Fourth order calculation conserve the total energy within 
the precision AE/E < 10~ 8 . 

Appendix C: Effects of windowing 

In this Appendix, we discuss reduction of white noise 
by windowing. For simplicity, here we identify time t and 
frequency uj with the rescaled variables t' and uj' . 

In general, the Fourier transform of a real-time func- 
tion fit) is defined by 

F{uj) = Jdt f{t)e wt (CI) 

From the definition, Parseval's identity 

/^l^)| 2 = jdt\f{t)\\ (C2) 

follows. In practice, we use f(t) at discrete time with 
time step At and upper bound i max - In this case, the 
integrals are replaced by discrete sums 

[ dt^At > ( C3 ) 

/At 

/s-*r- £ • <«> 

0<m<t max /At 

Assuming white noise, we consider the case |/ no ise(i)| 2 = 
5 2 at any discrete time and |F no i se (w)| 2 is independent of 
frequency uj. Then Parseval's identity implies 

\F noisc (u)\ 2 = S 2 t m , x At. (C5) 

Namely, the white noise in the frequency space increases 
as the sampling time interval t max is increased. The white 
noise could be reduced simply by taking smaller interval 
tmax, but it enhances spectral leakage. A better alterna- 
tive is to use the window function with smaller width, ef- 
fectively reducing the sampling time interval. In this pa- 
per, we use the Gaussian window function. The Fourier 



transform with the windowing is defined as 

F(w) = J dt j(t)e- e l^e wt . (C6) 
The couterpart of eq. (|C2I) becomes 

/ £'^ )|2 = J dt\f{t)?e- tV ° 2 ■ (C7) 
For the white noise with |/ no ise(i)| 2 = S 2 , and we find 

l^noiscHI 2 = 5 2 ^aAt. (C8) 
In our study, /noise (i) corresponds to the noise in the 
dynamical correlation function (Lij(f)L^(O)). sampling 
time step is^ At ~ 0.1, and S ~ 0.02. For a = 100, 
this gives \F noiBC (uj)\ ~ 0.2, which is about 2 % of the 
observed signal. In other words, this signal-to-noise ratio 
is achieved with a = 100. 

On the other hand, there is a side-effect of the win- 
dowing. In the absence of the staggered field h! = 0, the 
exact spectrum is proportional to S(uj — H), according to 
the general result in the absence of anisotropy. However, 
in Fig. [5J the paramagnetic resonance appears to have a 
finite width and the peak seems to be shifted to higher 
frequency uj' f=s 0.015 from the expected uj' = H' = 0.005. 

This is actually due to the the artificial linewidth cr _1 
introduced through the window function (|64[) . It actu- 
ally dominates the linewidth of the paramagnetic peak in 
Fig. [51 Once a finite width is induced, the peak in the dy- 
namical susceptibility is distorted by the factor 1 — e~^l T 
(which is proportional to uj at low frequency and give 
more weights to higher frequencies) in the definition, 
eq. (150)) . For a uniform field H' sufficiently larger than 
the resolution, we find the (resolution-limited) paramag- 
netic resonance peak at uj' ~ H', as expected. In Fig. [51 
we have shown the data for small H' so that the addi- 
tional broad peak can be observable in a wide range of 
frequency within the validity of the classical approxima- 
tion. The low- frequency peak in Fig. [5J is smoothly con- 
nected to the paramagnetic resonance peak at uj' ~ if', 
when the uniform magnetic field H' is increased; thus it 
can be identified with the paramagnetic resonance peak 
despite the apparent width and shift due to the frequency 
resolution in the calculation. In fact, the paramagnetic 
peak itself can be observed with a larger value of a. 
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